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Abstract 

The  velocity  data  from  a  direct  numerical  simu¬ 
lation  of  low  Reynolds  number  turbulence  in  an  open 
channel  have  been  used  to  compute  the  terms  in  the 
budget  equations  for  the  turbulence  kinetic  energy, 
the  dissipation  of  turbulence  kinetic  energy  and  the 
Reynolds  stresses.  The  budget  data  show  that  the  dis¬ 
sipation  rates  of  the  horizontal  components  of  the  tur¬ 
bulence  are  reduced  near  the  surface  while  the  dissipa¬ 
tion  of  the  vertical  component  remains  approximately 
constant.  The  data  also  show  that  the  pressure-strain 
term  is  the  dominant  producing  term  for  the  spanwise 
component  of  energy  in  the  near  surface  region  A 
model  for  this  behavior  valid  for  flows  exhibiting  ho¬ 
mogeneity  in  the  spanwise  and  streamwise  directions 
is  proposed  and  tested  against  the  data  In  general 
the  model  is  found  to  work  well  but  wider  testing  is 
necessary. 


n 

n. 

n^,  etc. 
Subacrtpts 


coordinate  direction 
xiUt/j/ 

1  -  Z-i 

Kronecker’s  delta 

isotropic  dissipation  function 

averaged  quantity 

kinematic  viscosity 

decomposed  pressure  strain,  Eq  (9i 

density 

shear  stress 

instantaneous  vorticity  vector 
instantaneous  vorticity  component 
terms  in  Lqs.  (6,7) 


1,2,3,  coordinate  directions 
value  at  free  surface 
value  at  wall 


Oij  =s  Reynolds  stress  anisotropy 

Aa,A,A'  ~  functions  of  a,j 

Oa.iai  etc.  =  terms  in  Eqs.  (5) 

Ca  —  model  constant  ¬ 
ly  =  gravitational  constant 

h  —  channel  height 

k  —  turbulence  kinetic  energy 

f*  =  v/ur,  viscous  length  scale 

p  =  fluctuating  pressure 

=  defined  in  Eqs  (8a, b) 

Rth  —  U,hli/,  Reynolds  number 

R"  =  Urh/u,  wall  Reynolds  number 

S,j  =  mean  strain  rate 

t'  =  u/xt^,  viscous  timescale 

U  =  instantaneous  velocity  vector 

U,  —  instantaneous  velocity  component 

=  fluctuating  velocity  component 
Ut  =  y/r^Tp,  friction  velocity 

u*  =  Ui/Ut 
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The  turbulent  flow  below  a  gas-liquid  interface 
plays  an  important  role  in  diverse  areas  ranging  from 
environmental  flows  and  industrial  mLxing  processes 
to  the  remote  sensing  of  ship  wakes.  The  near¬ 
boundary  influences  upon  transfer  and  diffusion  at  the 
interface  are  of  primary  concern  in  environmental  and 
industrial  applications,  whereas  remote  sensing  issues 
ultimately  involve  any  surface  motions  that  may  be 
detectable.  For  example,  the  two  most  common  and 
persistent  features  seen  in  synthetic-aperture-radar 
images  of  ship  wakes  are  bright  “narrow  vees”  and 
long  dark  “scars”,  which  may  be  a  result  of  surface 
Bragg  wave  generation  or  modification  through  inter¬ 
actions  with  near  surface  turbulence.  Common  to  all 
of  these  problems  is  the  need  for  a  better  understand¬ 
ing  of  the  structure  of  turbulence  below  a  free  surface 

For  several  decades  it  has  been  realized  that  the 
presence  of  a  free  surface  influences  the  evolution  of 
mean  velocity  and  turbulence  but  the  mechanisms 
have  not  been  completely  described.  Early  observa¬ 
tions  of  .Nikuradse*  showed  the  flow  in  straight  open 
channels  to  be  three-dimensional  and  that  the  max¬ 
imum  of  the  streamwise  mean  velocity  occurs  below 
rather  than  coincident  with  the  free  surface  Mor** 
recently,  the  studies  of  Ueda*  and  Komori  et  al  ^  for 
open  channel  flows  showed  that  the  eddy  viscosity  is 
significantly  attenuated  by  the  presence  of  the  free 
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surface.  lo  the  latter  paper  it  is  also  shown  that  in 
a  region  near  the  surface  the  surface-normal  ve¬ 
locity  fluctuations  are  diminished  while  the  fluctua¬ 
tions  in  the  plane  of  the  surfttce  are  increased,  he 
largest  increase  in  the  near-surface  region  is  in  the 
spaowise  component.  This  paper  also  indicates  that 
the  viscous  dissipation,  e,  has  vanishingly  small  nor¬ 
mal  gradient  near  the  free  surface.  A  similar  redistri¬ 
bution  of  the  turbulence  intensities  was  observed  by 
Thomas  and  Hancock.^  fn  their  work  a  moving  wall 
experiment  was  devised  such  that  the  wall  moved  at 
the  velocity  of  the  adjacent  turbulent  fluid  so  that  no 
velocity  gradients  and  shear  stresses  were  present  at 
the  wall,  conditions  similar  to  those  at  an  uncontami¬ 
nated,  wavelesa  free  surface.  Damping  of  the  velocity 
fluctuations  in  the  wall-normal  direction  accompanied 
by  an  increase  in  the  streamwise  fluctuations  was  ob¬ 
served.  The  spanwise  fluctuations  were  only  slightly 
increased.  In  recent  experimental  studies  by  Ramberg 
et  al.®  and  Swean  et  al.®'^  single-point  hot-film  mea¬ 
surements  of  the  velocity  correlations  were  obtained 
near  the  free  surface  in  a  jet  flow  These  measure¬ 
ments  also  showed  the  existence  of  a  thin  layer  near 
the  free  surfsuie  wherein  the  redistribution  of  turbu¬ 
lence  energy  occurred  rapidly  with  most  of  the  verti¬ 
cal  component  transferring  into  the  spanwise  compo¬ 
nent,  an  observation  similar  to  that  of  Komori  et  al. 
referenced  above. 

All  of  the  above  experiments  experienced  prob¬ 
lems  in  acquiring  data  very  near  the  surface  The  hot- 
film  studies  suffer  from  the  effects  of  probe  contam¬ 
ination  and  blockage  brought  about  by  the  intrusive 
sensor  near  the  boundary.  Laser  doppler  velocimetry 
methods  have  problems  due  to  reflection  and  refrac¬ 
tion  at  the  free  surface.  As  a  result  of  the  experimen¬ 
tal  problems  associated  with  obtaining  reliable  mea¬ 
surements  near  an  air-water  interface,  the  situation  is 
that  considerably  leas  is  known  about  the  character¬ 
istics  of  turbulence  near  a  free  surface  as  opposed  to 
flow  near  to  solid  walls. 

Due  to  the  interest  and  relevsmce  of  the  prob¬ 
lem.  a  direct  numerical  simulation  of  turbulent  open 
channel  flow  has  been  performed  by  Leighton  et  al  ® 
This  paper  is  an  analysis  of  time-averaged  data  from 
Leighton’s  calculation  with  the  goal  of  evaluating  and 
improving  turbulence  models  for  use  in  practical  cal¬ 
culations  of  the  near-surface  flow.  In  the  next  section 
a  brief  description  of  the  calculation  is  given 

2.  Direct  Numerical  Simulation 

The  incompressible  3D  Navier-Stokes  equations 
were  solved  for  initial  and  boundary  conditions  ap¬ 
proximating  a  turbulent  open-channel  flow  of  wa¬ 
ter  at  Rept  =  2340  based  on  the  channel  depth,  h. 
and  the  mean  steady  velocity  at  the  free  surface. 


U,  =  =  h).  The  governing  equations  were  re¬ 

cast  in  the  manner  suggested  by  Orzag  and  Patera^ 
and  implemented  by  Kim,  Moin  and  Moeer^°  for 
closed  channel  flow.  The  final  equation  system,  in 
which  the  pressure  has  been  eliminated,  consists  of  a 
4***  order  equation  for  the  vertical  velocity, 


l(Uxn)j 


\  dt  {dxl^dxir 

_J_  f  d 

dxi 

and  a  2"“^  order  equation  for  the  vertical  vorticity, 


dxs 


(u  X  n)i 


-  A(u  X  n)„  (2, 

where  all  variables  are  non-dimensionalized  by  h  and 
U,  and  bold-face  type  indicates  vector  quantities  with 
n  =  (  V  X  U).  Following  solution  of  Eqs  (1,2). 
the  streamwise  and  spanwise  velocity  components 
(L^,  L’a)  are  recovered  from  the  incompressibility  con¬ 
dition  and  the  definition  of  vorticity. 

The  equations  are  solved  after  they  are  Fourier 
transformed  in  the  streamwise  (zi)  and  spanwise  (za) 
directions  and  Chebyshev  transformed  in  the  vertical 
direction  (zj).  The  calculations  were  performed  on 
a  48  X  65  X  64  grid  in  zj.zj.za  respectively,  which 
allows  the  resolution  of  ail  essential  turbulent  scales 
without  resort  to  subgrid  models.  With  the  geome¬ 
try  scaled  by  the  channel  height,  the  vertical,  stream- 
wise  and  transverse  dimensions  of  the  channel  are  1 
4t  and  3)r/2.  respectively.  In  wall  units  the  domain 
IS  134f  X  1684f  X  632f,  where  f  =  h/R'  with 
R'  —  134  For  comparison  purposes  a  companion  cal¬ 
culation  was  performed  for  a  closed  channel  flow  For 
reasons  of  economy  this  calculation  was  at  half  the 
wall-normal  resolution  of  the  open  channel  case  and 
was  at  a  lower  Reynolds  number,  R"  =  125  Nev¬ 
ertheless  the  qualitative  behavior  of  the  data  for  all 
aspects  examined  was  identical  to  that  reported  by 
Mansour.  Moin  and  Kim  (henceforth  MKM).'*  As 
such  these  closed  channel  data  can  serve  in  certtun 
instances  to  compare  qualitatively  the  different  be¬ 
havior  in  the  open  and  closed  channel  simulations 


The  boundary  conditions  are  periodic  on  all  de¬ 
pendent  variables  in  the  streamwise  and  spanwise  di¬ 
rections.  No  slip  conditions  are  used  at  the  chan¬ 
nel  bottom  while  the  free  surface  is  approximated  as 
a  rigid  free  slip  surface  with  vanishing  shear.  The 
shear-free  rigid  lid  condition  is  an  approximation  to 
the  exact  free  surface  condition  which  is  valid  at  low 
Froude  number  for  a  surface  free  of  any 

contaminants  Leighton  et  al.*  have  estimated  the 
surface  displacement  s  psstenisn  from  the  resulu  of 
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10®  10'  10^ 
Fig.  1  Mean  velocity  profile  across  the  channel 


Fig  2  Profiles  of  turbulence  kinetic  energy  and 
Reynolds  stress  across  the  channel. 


the  simulation  and  using  the  channel  height,  h  = 
0  04m,  from  the  experiments  of  Komon  et  al  ^  For 
these  conditions  the  rms  surface  deflections  could  be 
expected  to  be  approximately  1.6  x  10"‘'m  (0  004/i) 
and  negligible  as  observed  in  the  laboratory  experi¬ 
ments.  For  later  reference  the  boundary  values  of  the 
dependent  and  derived  variables  at  the  wall  =  0) 
and  free  surface  {x^  =  1)  are: 

=t/,  =  f.r3  =  n,=  ^  =  0;X2  =  0,  (3) 

0X2 


and 


dUi 

dx2 


dUj  _ 
dzj 


dQj  _ 
dx2  ~ 


d^U2 

dx\ 


=  0;  Xt  =  1. 


(4) 


The  derivative  conditions  on  U2  arise  from  continuity 
considerations  at  the  respective  boundaries 

The  computer  code  used  in  the  simulation  was 
designed  and  developed  to  run  on  the  CRAY  X- 
MP/24  at  the  Naval  Research  Laboratory  Approxi¬ 
mately  10~®  seconds  per  timestep  per  grid  point  were 
required  for  the  simulation.  After  the  wall  shear  stress 
achieved  a  statistically  steady  behavior.  42  realiza¬ 
tions  of  the  instantaneous  vclcx:ity  data  were  saved 
during  a  time  interval  of  approximately  4000<'  where 
t*  =  u/u^.  Statistics  were  obtained  by  averagmg  in 
the  streamwise  and  spanwise  directions  and  in  time 

3.  Discussion  of  Results 


wall  laws  u+  =  x+  and  u+  =  2.5  In  a;  J  -I-  5.5  A 
best  fit  of  the  present  data  for  the  logrithmic  region 
is  u"*"  =  2.41nr2  "F  lower  slope  is  consis¬ 

tent  with  the  value  of  2. 43  found  by  Nezu  and  Rodi'^ 
for  open  channel  flows  over  a  Reynolds  number  range 
439  <  /?*  <  6139.  The  intercept  is  near  the  up¬ 
per  bound  (5.29  ±  0.47)  found  in  their  experiments 
and  is  probably  a  low  Reynolds  number  effect  ,"  Th** 
notable  difference  between  the  velocity  data  in  Fig  1 
and  closed  channel  behavior  is  the  absence  of  a  clearly 
defined  wake  region  in  the  outer  flow,  rather  the  Ion, 
law  IS  maintained  until  very  close  to  the  free  surface 
when  the  velocity  adjusts  to  the  vanishing  gradient 
condition. 

Figure  2  shows  normalized  turbulence  kinetic 
energy  k  =  -)-  Uj  u^)  and  the  Reynolds  shear 

stress,  U1U2,  which  has  been  further  scaled  for  plot¬ 
ting  by  a  factor  as  shown.  According  to  the  boundarv 
conditions  given  by  Eqs.  (3,4)  the  velocity  compo 
nents  and  pressure  can  be  expanded  about  the  free 
surface  as. 


«i  =  <Ji  -t-ciy2 

=  ijy  +  djy®  +  0(i/®) 

W3  =  as  +  C3J^  +  0(]^)  ' 

P  =  ap  +  cpv2+0(y®) 

where  y  is  defined  with  the  origin  at  the  free  surface 
Use_of  tjiese  expansions  and  averaging  results  m  F  = 
5(0^  -t-  a^)  +  O(y^)  near  the  surface  and  dk/dy  - 
at  the  surface  which  is  evident  in  the  figure. 


The  mean  velocity  normalized  by  Ur  is  shown 
in  Fig.  1,  Also  shown  by  the  dashed  lines  are  the 
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Fig  3  Profiles  of  aormaJ  stresses  across  the  channel 
and  comparison  with  data  of  Komori  et  al  ^ 

Figure  3  shows  the  three  components  of  k  from 
the  simulation.  The  symbols  are  the  experimen¬ 
tal  results  of  Komori  et  al.^  which  were  taken  at 
Reii  ss  3100.  There  is  qualitative  agreement  between 
the  experimental  and  numerical  results  Both  exhibit 
an  increase  in  the  horizontal  components  as  the  free 
surface  is  approached  and  the  increase  is  greatest  for 
the  spanwise  component.  The  computed  local  mini¬ 
mum  of  the  spanwise  component  occurs  further  from 
the  free  surface  than  in  the  streamwise  component 
which  is  also  consistent  with  the  experiments  This 
behavior  is  more  easily  recognizable  in  Fig  4  which 
contains  the  distributions  of  the  three  components  of 
the  turbulence  kinetic  energy  made  non-dimensional 
by  the  local  value  of  k.  Also  shown  in  this  figure 
are  the  data  from  the  closed  channel  calculation  ref¬ 
erenced  in  the  previous  section  and,  as  noted,  these 
data  are  presented  for  qualitative  comparison  only 
For  the  closed  channel  calculation  the  xj  =  1  bound¬ 
ary  corresponds  to  the  channel  centerline  Figure  4 
shows  that  near  the  free  surface,  xj  >  0  7,  most  of 
the  energy  from  the  vertical  component  is  transferred 
to  the  spanwise  component  with  only  a  small  increase 
in  the  horizontal  component.  This  contrasts  with  the 
behaviors  of  the  various  components  of  turbulence  en¬ 
ergy  in  the  closed  channel  simulation  where  the  rel¬ 
ative  interchange  of  energy  appears  to  be  primarily 
from  the  streamwise  component  to  the  vertical  in  the 
region  near  the  channel  centerline.  In  order  to  un¬ 
derstand  this  behavior,  the  budget  equations  of  the 
individual  velocity  correlations  have  been  examined 


Fig.  4  Profiles  of  normal  stresses  normalised  by  the 
local  value  of  turbulence  kinetic  energy.  Dashed  lines 
are  from  a  closed  channel  simulation. 

For  the  flow  under  consideration,  which  is  statis¬ 
tically  steady  and  homogeneous  in  the  spanwise  and 
streamwise  directions,  the  transport  equations  for  the 
one-point  velocity  correlations  are, 

=  0  =  +  n.;  +  r.;  +  +  D.,  -  .  (6) 

The  symbols  on  the  right-hand  side  of  F,q.  (6)  denote 
the  rates  of  production,  pressure  diffusion,  turbulent 
transport,  pressure  strain,  viscous  diffusion  and  dis¬ 
sipation.  respectively  The  explicit  representations  of 
these  terms  are: 


The  equation  obtained  by  taking  half  the  trace  of 
Eq  (6)  is  the  equation  for  the  turbulence  kinetic  en 
ergy,  k.  The  equation  for  the  trace  of  the  dissipation 
rate  tensor,  c  =  (<n+t5S+cas)/2,  is  given  by  Hanjalic 
and  Launder*^  as, 
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Fig.  5  Distribution  of  the  terms  in  the  budget  of 
in  the  upper  half  of  the  channel. 


Fig.  6  Distribution  of  the  terms  in  the  budget  of 
in  the  upper  half  of  the  channel. 
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The  first  four  terms  are  production  terms  (P/  to  P*). 
while  terms  five  through  eight  are  pressure  transport 
(n,),  turbulent  transport  (T<),  viscous  diffusion  (D,) 
and  dissipation  (Y),  respectively. 


The  terms  in  the  budget  equations  for  the  three 
normal  stresses  and  the  dissipation  rate  are  shown  in 
Figs,  S-8.  All  terms  in  Eqs,  (6,7)  have  been  normal¬ 
ized  by  and  the  budgets  are  displayed  only  for 
the  upper  half  of  the  channel  nearest  the  free  surface 
Since  the  production  of  turbulence  is  much  lower  than 
in  the  high  shear  region  close  to  the  solid  wall,  the  in¬ 
dividual  terms  in  the  equations  are  typically  an  order 
of  magnitude  lower  than  their  corresponding  values 
near  the  wall.  The  near- wall  data  are  available  in 
Leighton  et  al.*  and  are  very  similar  to  the  data  of 
MKM  " 


Figure  5  shows  that  away  from  the  surface  all 
terms  in  the  budget  have  the  same  relative  im¬ 
portance  except  Dll,  the  viscous  diffusion  Moving 
toward  the  surface,  the  production  rate  vanishes  with 
the  mean  velocity  gradient.  At  the  wall,  the  viscous 


terms  balance  the  turbulence  transport  and  the  pres¬ 
sure  strain,  the  latter  having  become  a  slight  positive 
contributor  to  the  budget  very  near  the  wall.  The 
ul  balance  in  Fig.  6  is  relatively  more  complex  than 
that  for  For  this  component  the  magnitude  of  the 
budget  terms  near  the  free  surface  are  only  reduced  In 
about  one-half  from  their  values  near  the  solid  wall 
Near  the  free  surface,  the  aaymtotic  behavior  of  thf- 
various  terms  can  be  determined  by  using  Eqs  (5i  as 

T22  =  +  •  ■ 

[122  =  “20p62  ~  6(Cp62  +Opd2)j/^  -t- 
•^22  =  2op62  -f  2(cpb2  -F  3opd2)y^  -F 
D22  —  262^2  -F  2462^21/^  -F  ■  ■  ■ 

<22  —  262^2  "F  12621^251^  -F  ■■ 


It  is  seen  that  at  the  free  surface  (22  balances  D22  and 
the  two  pressure-velocity  terms  cancel.  Note  that  the 
pressure  strain  has  rapidly  become  a  consuming  term 
in  the  near-surface  region  whereas  it  had  been  a  ma 
jor  producer  in  the  budget  equation  in  the  outer  flow 
This  is_in  contrast  to  the  behavior  shown  in  Fig  7 
for  the  component.  In  this  case  *33  increases  near 
the  surface  and  at  the  surface  is  considerably  mor*- 
of  a  source  for  than  is  #11  in  the  uf  budget  This 
largely  explains  why  the  transverse  component  of  tur 
bulence  kinetic  in  Fig,  4  is  increased  relatively  mor» 
so  than  the  streamwise  component  as  Oj  approaches 
zero  at  the  free  surface.  In  the  next  section  a  mole; 
for  this  behavior  is  proposed. 
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Fig  7  Distribution  of  the  terms  in  the  budget  of  ti| 
in  the  upper  half  of  the  channel. 

A  curious  feature  contained  in  Figs  5  and  7  is 
the  behavior  of  «,*  as  the  surface  is  approached  Mov¬ 
ing  toward  the  free  surface  the  dissipation  rates  ex¬ 
hibit  a  sharp  drop  in  magnitude  in  the  upper  5-10%  of 
the  channel.  The  <32  on  the  other  hand  shows  a  very- 
slight  increase  in  magnitude  and  could  well  be  approx¬ 
imated  as  constant  in  this  region  This  behavior  is 
contraiy  to  standard  modeling  assumptions  near  the 
free  surface.  Hossain  and  Rodi*^  and  later  Naot  and 
Rodi‘®  have  assumed  that  in  most  respects  other  than 
the  vanishing  of  the  surface-normal  velocity  compo¬ 
nent,  the  free  surface  behaves  like  a  symmetry  plane 
The  exception  is  the  presumed  behavior  of  the  dis¬ 
sipation  rate  for  turbulence  kinetic  energy  which  is 
expected  to  increase  near  the  surface.  This  is  based 
on  the  assumption  that  the  macro-length  scale  of  the 
turbulence  (£  oc  fc*  */e)  is  reduced  by  the  presence  of 
the  boundary.  This  scale  does  not  become  zero  since 
it  reflects  the  fluctuating  motion  in  all  three  direc¬ 
tions  and  the  horizontal  extent  of  the  eddies  is  not 
restricted. 

The  terms  in  the  budget  equation  for  the  dissi¬ 
pation  rate  are  shown  in  Figure  8.  In  the  upper  por¬ 
tion  of  the  channel  the  first  three  production  terms 
in  Eq.  (7)  are  small  and  have  been  lumped  together 
as  shown  Until  very  near  the  free  surface  the  pro¬ 
duction  by  turbulence.  P*.  largely  balances  the  dis¬ 
sipation  term.  Very  near  the  surface  the  dominant 
terms  are  the  viscous  diffusion  and  the  dissipation, 
each  exhibiting  very  large  gradients  of  opposite  sign 
near  the  boundary  Reconsidering  Figs  5  and  7.  it 
IS  seen  that  in  these  cases  also  the  rapid  variation  in 


Fig.  8  Distribution  of  the  terms  in  the  budget  of  the 
turbulence  kinetic  energy  dissipation  rate  in  the  upper 
half  of  the  channel 

the  dissipation  term  appears  to  be  at  least  partially 
offset  by  the  variation  in  the  diffusion  term.  It  should 
be  remembered  that  f„  (or  (}  is  not  the  actual  dis¬ 
sipation  of  turbulent  energy  for  inhomogeneous  flows 
although  it  does  approximate  the  total  dissipation  for 
high  Reynolds  number.  The  particular  terms  D„  and 
€„  have  arisen  from  the  combination  of  the  actual 
dissipation  rate  with  the  rate  of  work  by  the  viscous 
shear  stresses  of  the  turbulence.^®  In  flows  far  from 
solid  walls  the  viscous  diffusion  is  generally  neglected 
and  as  such  the  modeled  dissipation  rate  implicitly 
mode''  the  work  term  Figure  9  shows  the  balance  of 
turbulence  kinetic  energy  obtained  from  the  trace  of 
Eq  (6)  The  viscous  terms  have  been  added  and  to¬ 
gether  they  balance  the  transport  terms  at  the  wall  It 
is  seen  that  the  total  viscous  term  varies  only  slightly 
near  the  boundary  and  might  be  easier  to  model  Fu¬ 
ture  work  will  re-process  the  velocity  data  to  deter¬ 
mine  the  actual  dissipation  term  separately  from  th'* 
work  term. 

•I  Pressure-Strain  Rate  Model 

The  above  data  show  that  the  pressure-strain 
term  is  a  key  contributor  to  the  redistribution  of  nor 
mal  stresses  near  the  boundaries.  Leighton  et  al  ® 
have  decomposed  the  fluctuating  pressure  into  a  slow 
pressure,  p''*,  a  'feist'  pressure.  p'^U  and  a  Stok^' 
pressure,  p'*',  in  the  manner  suggested  by  MKM  ' 
These  components  satisfy  the  equations. 
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Fig.  9  Distribution  of  the  terms  in  the  budget  of  tur¬ 
bulence  kinetic  energy  m  the  upper  half  of  the  chan¬ 
nel. 
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The  pressure  data  from  the  solution  of  Eqs  (8a-b) 
have  been  used  to  decompose  the  pressure-strain  term 
into 


Figure  10  shows  this  splitting  in  the  upper  portion  of 
the  channel  for  each  of  the  components  The  Stokes 
term  is  only  of  consequence  very  near  the  solid  wall 
and  IS  not  plotted.  It  is  seen  that  the  slow  or  return 
terms  are  dominant  in  the  region  plotted  Only  in 
the  <^33  component  is  the  fast  term  the  larger  of  the 
two  and  for  the  <^22  term  the  fast  component  is  nearly 
zero  over  the  whole  domain  shown. 


Fig  10  Distribution  of  the  rapid  and  return  contri¬ 
butions  to  the  pressure-strain  rate  for  each  diagonal 
component 


The  distributions  of  shown  in  Fig  10  differ 
markedly  at  the  free  surface  relative  to  their  behav¬ 
ior  at  a  solid  wall.^*  Most  importantly  they  do  not 
vanish  in  any  component.  Clearly  some  sort  of  sur¬ 
face  proximity  effect  is  needed  in  order  to  model  the 
rapid  variations  exhibited.  Generally  the  total  term 
is  modeled. 
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where  is  some  variant  of  Rotta’s''  model. 


♦Ij' =  -Gea.;,  (101 

and  the  rapid  term  includes  at  least  the  isotropization 
of  production  term.  (Pij  —  ^PtS,j).  Typical  of  these 
models  and  one  that  is  borrowed  from  in  the  rurreni 
study,  is  that  due  to  Launder.  Reece  and  RodiO^ 
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where. 


„  _ dUt  _ dUt 


and  C2  =  0  4.  The  wall  terms  are  usually  modeled 
in  a  form  analagous  to  Eqs.  (10,11)  but  with  different 
coefBcients  and  including  a  damping  function  The 
Launder.  Reece  and  Rodi  form  is. 
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the  horizontal  componenU  will  tend  to  approach  this 
state  at  a  free  surface.  The  available  experimental 
data  as  well  as  the  current  simulation  data  support 
this  conjecture.  If  this  is  the  case  then  the  vanish¬ 
ing  of  A3  can  be  used  as  a  detector  for  free  surface 
proximity  effects. 

A  model  based  on  Lumley’s  suggestion  and  used 
as  the  basis  of  the  current  model  is  that  due  to  Laun¬ 
der  and  Shima.^°  In  their  model  is  given  by 
Eq  (10)  with  the  coefficient, 

Cl  =  2.58AA,“{1  -exp((0.0067*Vi^O^)}  n,3i 

The  rapid  term  is  given  by  the  first  term  in  Eq  n.  1 1 
with  the  coefficient,  SfCj  +  S)/!!,  replaced  with  (  ’2  = 
0  T.o.q  ^  The  Launder  and  Shima  model  also  contain' 
additional  explicit  expressions  for  the  wall  terms  m 
Eq  (9)  that  are  not  repeated  here  since  they  ar>>  nm 
used 


Fig  1 1  Distribution  of  the  second  and  third  invariants 
of  the  Reynolds  stress  anisotropy  tensor  across  the 
channel. 

This  particular  damping  function  is  singular  at  a  free 
surface  since  k  does  not  vanish.  A  virtual  origin  could 
be  used  and  has  been  employed  by  Naot  and  Rodi** 
An  alternative  is  to  use  a  term  based  on  the  surface- 
normal  velocity  such  as  (u^)  ^k/dn.  which  does 
not  become  singular  and  is  more  appealing  on  phys¬ 
ical  grounds.  However,  the  approach  taken  in  this 
study  is  to  make  use  of  the  properties  of  the  sec¬ 
ond  and  third  invariants,  Aj  and  Aa  of  the  stress 
anisotropy  tensor.  a,f. 

A2=a,,(3,;  ;  A3  =  a,;a;tat,  (12) 

.\s  pointed  out  by  Lumley’®.  if  one  component  of 
velocity  vanishes,  then  the  difference,  Aj  -  .43,  be¬ 
comes  the  constant  value  8/9  irrespective  of  the  be¬ 
havior  of  the  other  two  velocity  components  In  this 
case  the  function  A  =  1  —  9(A2  —  A3)/8  =  0  in  the 
regions  where  the  turbulence  becomes  locally  two- 
dimensional.  At  the  free  surface  the  pressure-strain 
rates  do  not  vanish  so  it  is  necessary  to  modify  this 
approach.  Figure  11  shows  the  variation  of  .-It  and 
.43  computed  from  the  simulation  data  Near  both 
boundaries  Aj  —  A3  — »  8/9  as  037  —  —2/3  Near  the 
free  surface  it  is  observed  that  A3  becomes  negative 
This  IS  very  close  to  the  17  location  where  033  becomes 
positive  (see  Fig.  4).  For  flows  in  which  the  non-zero 
velocity  components  remain  uncorrelated  as  the  third 
component  vanishes.  t..e  vanishing  of  A3  corresponds 
exactly  to  the  vanishing  of  one  of  the  remaining  re¬ 
maining  a„  Equation  (12)  shows  that  for  such  a  flow 
.43  reaches  a  minimum  of  -2/9  when  an  =  033  It  is 
plausible  to  assume  that  m  the  absence  of  boundaries 


In  the  current  work  an  effort  has  been  made 
to  construct  a  model  with  a  minimum  of  parame¬ 
ters  In  Eq  (13).  both  the  coefficient.  A.  and  the 
term  in  braces  containing  the  Reynolds  number  act 
to  damp  the  pressure-strain  rate  near  a  solid  wall 
Equations  (8a,b)  are  only  dependent  upon  viscosity 
thru  the  boundary  condition  on  the  equation  for  p'* 
and  .MKM  have  shown  that  for  low  Reynolds  flows  of 
the  type  under  study  here.  is  small  even  near  the 
wall  The  term  in  braces  has  thus  been  omitted  in 
the  current  study  Exploratory  calculations  retaining 
the  term,  and  using  the  statistics  from  the  simulati'-n 
have  shown  that  the  term  has  little  effect  on  mod-^l 
performance  and  in  fact  the  pressure-strain  rate  data 
are  better  correlated  with  its  omission  Although  th** 
concern  here  has  not  been  model  performance  ne.a.'' 
the  solid  wall,  this  study  has  also  retained  the  third 
term  in  Eq  (11)  since  its  inclusion  substantially  im 
proves  the  behavior  in  the  ^77  component  near  the 
wall  The  basic  form  of  the  current  model  then  he 
comes. 


<!>;)'  =  -C,Af a,;  ;  Cl  =  1.3A2^V  C-l 


(Q^8) 

ri~' 


-I  P.5,, ) 


in 


The  numerical  value  in  the  definition  of  C]  above  has 
been  adjusted  from  its  former  value  to  account  for 
the  factor  of  approximately  1/2  that  the  retention  of 
the  exponential  term  of  the  original  model  would  havp 
contributed  in  the  outer  flow.  The  form  adopted  f  t 
the  free  surface  region  is, 


-I-:;  .  =  -Ci(A  -  A')e(2.2a,;  -I- 9.8(a,*ai.  -  A2/3)i 
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F'g  12  Comparison  of  the  modeled  pressure -strain 
rate  with  the  data  from  the  simulation 

where  C]  is  given  by  Eq.  (14)  and  .4'  =  1  -  9/8(.43  - 
I.43I)  As  discussed  earlier  this  term  will  only  be¬ 
come  active  when  the  anisotropy  in  one  of  the  hori¬ 
zontal  components  vanishes  while  the  vertical  veloc¬ 
ity  18  damped  near  the  free  surface.  The  constants 
2,2  and  9  8  have  been  determined  so  as  match  the 
overall  level  of  the  data  m  the  surface  region  The 
2"^  order  term  has  been  shown  by  Speziale  et  alT' 
to  occur  in  a  formal  expansion  for  <l>,;  and  it  is  used 
in  the  current  formulation  to  achieve  more  closely  the 
separation  in  ^^33  ttnd  observed  in  the  simulation. 
Near  the  free  surface  the  anisotropy  in  the  U3  com¬ 
ponent  of  velocity  is  relatively  small  compared  to  au 
so  10  first  order  in  anisotropy  the  model  predicts  that 
*^11  >  ^33  The  second  order  term  is  small  in  both 
components  but  has  a  much  greater  relative  effect  on 
the  «3  component. 

Figure  12  shows  the  model  performance  when 
computed  with  the  simulation  data  which  are  rep¬ 
resented  by  the  symbols.  The  lines  are  computations 
with  Eqs  (14-16).  The  incorrect  near-surface  asym- 
totic  behavior  is  mostly  due  to  the  sharp  variation 
in  (  discussed  in  the  previous  section  Incorporat¬ 
ing  the  proper  asymtotic  behavior  near  the  free  sur¬ 
face  will  be  a  topic  for  further  research  Near  the 
solid  wall  the  current  model  performs  reasonably  well 
and  markedly  better  than  the  original  Launder.  Reece 
and  Rodi  model  which  was  analyzed  by  MKM  "  The 
asymtotic  behavior  near  the  wall  is  not  correct  and 
this  must  be  examined  as  well.  Without  the  inclusion 
of  the  term  containing  B,j  in  Eq.  (15)  the  change  in 
sign  in  the  near  wall  region  for  ♦33  is  not  achieved  but 


rather  the  modeled  distribution  is  uniformly  positive 
until  vanishing  at  the  wall. 

5.  Concluding  Remarks 

The  simulation  data  show  that  there  is  a  prefer¬ 
ential  redistribution  of  turbulence  energy  to  the  span- 
wise  component  of  energy  as  the  normal  component 
IS  damped  at  the  free  surface.  This  is  in  accordance 
with  the  still-limited  experimental  observations.  The 
budget  equations  show  that  the  pressure-strain  rate 
term,  particularly  the  return  term,  is  a  key  contrib¬ 
utor  to  this  behavior.  The  budgets  also  show  that 
the  isotropic  part  of  the  dissipation  of  turbulence  en- 
ery  decreases  rapidly  very  near  the  free  surface  which 
18  contrary  to  current  modeling  asumptions  which 
assume  it  to  increase  to  account  for  reduced  levels 
of  eddy  viscosity  observed  near  the  free  surface  .\ 
model  for  the  near-surface  pressure-strain  term  has 
been  proposed  and  shown  to  correlate  the  .simulation 
data  fairly  well.  In  its  current  form  the  model  is  lim¬ 
ited  to  flows  in  which  the  correlation  of  the  horizontal 
velocity  components  is  small  while  the  vertical  com¬ 
ponent  vanishes.  Much  wider  testing  is  necessary  to 
determine  when  this  condition  exists.  Future  efforts 
will  consider  this  point  along  with  questions  raised  re¬ 
garding  the  dissipation  of  turbulence  at  the  surface 
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